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Abstract 



Scaling properties of an interface representation of the critical contact process 
are studied in dimensions 1-3. Simulations confirm the scaling relation 
Pw = 1 — Q between the interface-width growth exponent flw and the exponent 
9 governing the decay of the order parameter. A scaling property of the height 
distribution, which serves as the basis for this relation, is also verified. The 
height-height correlation function shows clear signs of anomalous scaling, in 
accord with Lopez' analysis [Phys. Rev. Lett. 83, 4594 (1999)], but no 
evidence of multiscaling. 

PACS numbers: 64.60.Ht, 05.40.-a, 05.10.Ln, 05.90.+m 
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I. INTRODUCTION 



Scaling and criticality in nonequilibrium systems continue to be of great interest in 
statistical physics. Among the various classes of systems that have been subject to intensive 
study are models of growing interfaces |l|-f|] , and absorbing-state phase transitions, typified 
by directed percolation (DP) |5|-§. Absorbing-state phase transitions have been linked 
to self-organized criticality (SOC) [PHTlj, as have driven interface models [[i2H T5||. The 



latter connection is established by defining a height variable hi(t) at each site of a sandpile 
model; at each time interval, the height increases by one unit at each active (toppling) 
site. It turns out that the interface of the Bak-Tang-Wiesenfeld model |TB[ is described 



by an Edwards- Wilkinson equation with columnar noise [14,15|. In view of the connections 



between absorbing-state phase transitions, SOC, and surface growth, it is of interest to study 
the dynamics of the interface representation of a simple model in the DP class. Precise results 
on the scaling properties of the DP interface representation should prove useful when trying 
to assign interface models (or height representations of other models, such as sandpiles) to 
universality classes. 

In this paper we examine the dynamics of the contact process in dimensions 1-3, 
studying the interface width and the height-height correlation function, as well as the height 
probability distribution. In Sec. II we define the model and then present a brief discussion of 
the associated continuum equation, and of a scaling theory. Numerical results are presented 
in Sec. III. Sec. IV contains a brief summary. 



II. MODEL 

The contact process (CP) is a simple particle system (lattice Markov process) exhibiting 



a phase transition to an absorbing (frozen) state at a critical value of the creation rate ||T7 . 
This model belongs to the universality class of directed percolation f5[ and Reggeon field 
theory JI8fl , and is pertinent to models of epidemics ]7|], catalysis and damage spreading 
f2(| , among many others. In the CP each site of the hypercubic lattice Z d is either vacant 
or occupied by a particle. Particles are created at vacant sites at rate Xn/2d, where n is the 
number of occupied nearest-neighbors, and are annihilated at unit rate, independent of the 
surrounding configuration. The order parameter is the particle density p; the vacuum state, 
p = 0, is absorbing. As A is increased beyond A c , there is a continuous phase transition from 
the vacuum to an active state; for A > A c , p ~ (A — AJ' 3 in the stationary state. In one 
dimension, A c ~ 3.297848. 

There are a number of ways (equivalent as regards scaling behavior), of implementing 
the CP in a simulation algorithm; this work follows the widely used practice of maintaining 
a list of all occupied sites. In this study the initial condition is always that of all sites 
occupied. Subsequent events involve selecting (at random) an occupied site x from the N p 
sites on the list, selecting a process: creation with probability p — A/(l + A), annihilation 
with probability 1 — p, and, in the case of creation, selecting one of the 2d nearest neighbors, 
y, of x. (The creation attempt succeeds if y is vacant). The time increment At associated 
with an event is 1/N P , where N p is the number of occupied sites immediately prior to the 
event. A trial ends when all the particles have vanished, or at the first event with time > t m , 
a preset maximum time. 
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The important scaling laws pertinent to the critical contact process on a lattice of L d 
sites, starting with all sites occupied are: (1) the mean survival time r oc Z/li/^, and (2) 
the average particle density decays as a power law, pit) oc t~ e for 1 < t < r. (In practice, 
the power law is found already for t ~ 2.) In one dimension, — 1.5808 and 9 ~ 0.1595 

Occupied sites represent activity, which spreads from site to site. (The absence of any 
activity corresponds to the absorbing state.) In a growing surface or driven interface, activity 
corresponds to the motion of the interface. We therefore define the height hi(t) at site i as the 
amount of time (up to time t) that site i has been occupied. In our numerical studies we use a 
real- valued height (recall that time is not restricted to integer values in our implementation). 
By keeping a record of the last time £j at which the state of site % changed, we are able to 
evaluate hi(t) at any moment in the simulation. (While the results reported here are for real 
h, we find the same scaling properties for integer h.) The surface hi(t) may be thought of 
as a driven interface. Since the critical contact process must eventually enter the absorbing 
state, the interface, in this analogy, will eventually be pinned. 

The large-scale properties of the CP and related models such as DP can be described via 
a field theory (so-called Reggeon field theory) framed in terms of a coarse-grained density 
p(x, t) > jTj?]. Retaining only relevant terms, the stochastic pde for p(x, t) reads 

^ = V 2 p-ap-bp 2 +v(x,t) . (1) 

Here ry(x, t) is a Gaussian noise with zero mean and autocorrelation 

77(x,t)77(x',t') = r P (x, t)5 d (x - x')5(i - . (2) 
In our continuum description, the height is given by 

h(x, t) = /"* dt'p(x, t') . (3) 
Jo 

Integrating Eq. ([j]) from time zero to time t, we obtain 



dh 
~dt 



V*h-ah-bjU^\ dt' + C(*,t) 



(4) 



with the noise autocorrelation 



C(x, t)C(x', f ) = T6 d (x - xO^x, *<) , (5) 

where t < = min(t, t'). Thus the equation governing the height includes a nonlinear memory 
term and a noise with nonvanishing correlations between different times. This equation 
does not seem to shed much light on the scaling properties of the interface. Rather, the 
relation between h and the density in the contact process, Eq. ([|), yields some properties 
that are not immediately obvious from Eq. (HI), for example, that dh/dt > 0, and that the 
nonlinear term is relevant for d < d c — 4. In any event, Eq. (||) does serve to point up 
some differences between the CP interface and conventional surface-growth models such as 



the Edwards- Wilkinson [22] or Kardar-Parisi- Zhang equations f23 |. First, the linear "drive" 
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term (— a/i(x, t) with a < in the active regime) is proportional to the local height and so 
cannot be transformed away. Second, in the active state, the moduli of the last three terms 
on the r.h.s. of Eq. grow without limit, suggesting that V 2 h does as well, so that the 
width of the active phase never saturates. 

We turn now to a simple scaling analysis of the interface. At its basis lies a scaling 
hypothesis for the probability density p(h; t) of the height h (at any lattice site) at time 
t: the time-dependence of this density enters only through the mean height h(t). The 
conjectured scaling property is 

p(h;t) = J-V(h/h(t)) , (6) 

where the scaling function V > with / V(u)du = 1, and the prefactor guarantees normal- 
ization. It follows that the moments of h all scale with the mean height, h n = Ti !T [h(t)}' n with 
vP- the n-th moment V. In particular, the mean-square width 

W 2 {t,L) = var(/i) ~ [h(t)} 2 , (7) 

if Eq. (H) holds. On the other hand, we have that for times t < r in the critical contact 
process, 

h(t) = f dt 1 'p(t') ~ t 1 - , (8) 
J o 

which immediately implies that W 2 ~ t 2/3w with f3y/ = 1 — 9. 

In surface growth studies the crossover time is expected to scale as t x ~ L z ; for a process 
in the DP universality class, we may write z = v\\fu±_ (clearly t x and r should scale with the 
same exponent). Then the roughness exponent a, defined via W 2 (t,L) = W 2 at (L) ~ L 2a 
for t ^> t x , is given by the scaling relation a = fiw z = (1 — ®) v \\I v l- It is perhaps worth 
stressing that the expressions relating a and fiw to DP exponents depend on the validity of 
the scaling hypothesis, Eq. ([]), which should be tested. Inserting the known DP exponent 
values, scaling theory yields a ~ 1.3287, 0.97, 0.51, and zero for dimensions 1, 2, 3, and 
4, respectively. In other words, the interface associated with the critical contact process is 
super- rough in one dimension, and asymptotically flat in d > d c = 4, where 6=1. 

The family of height-difference correlation functions 

G g (r,t) = \h(x,t) - h(x + r,t)\i , (9) 

are also much studied in surface growth processes. Starting with a flat interface at t = 0, 
we expect power- law growth, G q ~ r gaq , for r < £(£) ~ t 1 ^, the time-dependent correlation 
length. If a q depends on q the interface is said to be multi-affine. For r > £, G q will saturate; 
in particular, G2 should approach W 2 (t,L) ~ t 2l3w ~ (, 2a for r ^ ^ and t < t x . Since £ 
is the only length-scale relevant to correlations at short times (i.e., for t < t x , so that the 
system size L does not come into play), it is reasonable to expect the scaling form 

G 2 (r,t)=e a Q(r/0, (10) 

where the scaling function Q{x) ~ x 2ot2 for small x and is constant for large x. The case 
a?2 = a is referred to as "conventional" scaling while a 2 < a is characterized as "anomalous" 
scaling 
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III. SIMULATION RESULTS 



We determined the square width W 2 (t,L) for rings of L = 500, 1000, 2000, and 5000 
sites, in samples of 2000, 1000, 1000, and 400 trials, respectively. The maximum time ranged 
from about 1.6 x 10 5 for L = 500 to 8.8 x 10 6 for L = 5000. All simulations were performed 
at the critical point, A c ~ 3.297848. In the contact process, interfacial properties can be 
studied over the full sample, or over a sample restricted to those trials that survive to time t. 
When a trial enters the absorbing state, W 2 naturally remains fixed for all subsequent times, 
and since all trials do eventually reach the absorbing state in the critical contact process, the 
interface width saturates for large t. If we restrict the sample to surviving trials, however, 
there is no saturation. It is important to note that the same scaling laws apply in either 
case. (The power-law growth regime, for example, corresponds to times < r, for which all 
trials still survive.) 

Fig. 1 shows a series of snapshots of the interface in a single trial with L = 200, at 
intervals of 5000 time units. The progressive roughening of the interface, without evidence 
of saturation, is apparent. Fig. 2 is a scaling plot of the square width, averaged over all trials, 
i.e., W 2 (L, t)/L 2a versus t = t/L z , using the exponents a = 1.328 and z = 1.5808 expected 
for DP in 1+1 dimensions. There is a near-perfect data collapse for t > 10~ 3 . The power-law 
portion of the graph (10~ 3 < t < 0.05) yields the growth exponent @w — 0.839(1), in good 
agreement with the value 1 — 9 = 0.8405 expected on the basis of the scaling argument. 
(Figures in parentheses denote statistical uncertainties.) Analysis of the saturation width 
yields W 2 at (L) ~ L 2a with a = 1.325(15). 

In Fig. 3 we test the scaling assumption for the height probability distribution p(h) by 
plotting t°- 84 p(h) versus h/t°' m . (Recall that h is expected to grow oc t 1 ^ = t 0,84 .) In this 
system of L = 1000 sites, there is a near-perfect data collapse, in accord with Eq. @, for 
times between 500 and 10 4 . For t = 2 x 10 we begin to note a departure from scaling, which 
becomes more pronounced at later times. Note that t = 2 x 10 4 corresponds to Int = — 1, 
which is where W 2 begins to depart noticeably from a power law in Fig. 2. Analysis of p(h) 
for L = 5000 yields similar results. The form of the scaling function for a given reduced time 
t is independent of the system size, as shown in Fig. 4, which compares height probability 
distributions for systems of 1000 and 5000 sites at times corresponding to the same t. 

The height-height correlation function G2 in a system of 1000 sites is shown in Fig. 5, 
which is a double- logarithmic plot of G\ = G2/t 2/3w versus r* = r/t 1 '". There is a perfect 
data collapse for t = 500,.. .,10 4 (reduced times 10 -2 < t < 0.2), with a power-law portion 
Q 2 r 2a ' 2 with «2 — 0.625. At later times does not collapse, principally because the 
square width has begun to saturate (it no longer grows oc t 2f3w ). The value of a 2 is insensitive 
to system size. For L = 5000 we obtain a 2 = 0.623(2) fort = 2xl0 5 , and 0.644 for t = 5xl0 5 . 
Thus we may, with a high degree of confidence, adopt the estimate a.2 = 0.63(3), clearly 
much smaller than the roughness exponent a = 1.33 found in the analysis of the square 
width, indicating that this one-dimensional system exhibits anomalous scaling. 

Recently, Lopez argued that anomalous surface roughening is associated with a diverging 
height gradient |24j| . In view of the growing spikes evident in the profiles shown in Fig. 1, 
this association seems very probable in the present instance. Indeed, the mean-square height 
gradient 

S (t)=ww. (11) 
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diverges as a power law (see Fig. 6). We find s(t) ~ t 2K , with re = 0.439(1) for L = 1000, 
0.4335(4) for L = 2000, and 0.4337(4) for L = 5000; we adopt re = 0.4336(4) as our best 
estimate. From Lopez' analysis, one expects 

a q = a — zk. (12) 

Inserting our result for re and the DP values a = 1.32867(14) and z = 1.5808(1) in the 
r.h.s., we obtain a q = 0.643(1), in good agreement with the result found from analysis of 
the correlation function. 

We also studied the generalized height-height correlation function G g , Eq. for q = 
1/2, 1, 3/2, 2 and 3, in a system with L = 1000, at t = 10 4 , and found that the functions 
[Gg] 1 / 9 are identical, and give a q ~ 0.62. We may therefore conclude that the DP interface 
is self-affine not multi-affine. 

In two dimensions we studied systems of up to 256x256 sites at the critical point, A c ~ 
1.6488. The curves for W 2 again show a good collapse (see Fig. 7), and the derived exponents 
are in very good agreement with the expected values. Specifically, we find a = 0.97(1), (3w = 
0.550(5), and z = 1.765(10). Scaling relations combined with known DP exponent values 
yield a = 0.970(5) and (3 = 0.549(2), while current best estimates give z = v\\/v± = 1.766(2) 
for DP in 2+1 dimensions [25-p8|. 

For smaller system sizes we observe a transient in W 2 and all other measured quantities, 
which does not appear in one dimension. Anomalous roughening is also seen in two and three 
dimensions. Fig. 8 shows s(t) growing ~ t 2K with re = 0.33(1). The height-height correlation 
function exhibits a good collapse, as shown in Fig. 9; the initial power-law growth yields 
a 2 = 0.385(5); thus the scaling relation Eq. flT2D is well satisfied: a 2 - a + zk = 0.00(3). 

We determined the interface exponents in three dimensions (A c ~ 1.3169, system size 
50 3 sites, maximum time 10 4 ), though to somewhat lower precision, owing to the larger 
computational demand. The scaling relation yields j3w = 1 — 5 = 0.274(1), while our 
simulation results for W 2 give j3w = 0.27(1). In three dimensions we find a = 0.51(1), 
z = 1.90(5), and re = 0.22(1). Together with Eq. ([12]), these yield a 2 = 0.09(2). Our results 
for critical exponents are collected in Table |TV|. 



IV. SUMMARY 

Defining an interface representation for the contact process by analogy with similar 
representations for sandpile models, we verified the expected scaling relation j3 w = 1 — 9 
in dimensions 1-3, and the scaling property of the height probability distribution in one 
dimension. The local roughness exponent, a 2 , is smaller than the global value, a, indicating 
anomalous surface growth. This anomalous scaling is attended by a diverging local slope, 
s(t) = (V/i) 2 ~ t 2n . Our results for re are consistent with the scaling relation Eq. (|i~2"D 
derived by Lopez. There is, on the other hand, no evidence of multi-affinity in this process. 

An interesting point is that the process continues to exhibit anomalous scaling for d = 2 
and 3, even though a < 1 in these cases. While it was pointed out some time ago that a > 1 
implies anomalous growth ||29|| , the latter appears to be an intrinsic feature of the contact 
process (and, by extension, of other models in the DP universality class) below d c . 

Finally, we note that the anomalous roughening analysis introduces two new critical 
exponents, a 2 and re, and only one new scaling relation between them. We therefore have a 
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new independent exponent, k for instance, that is not related to the standard DP exponents 
in any way. A very interesting theoretical task is that of computing k in an epsilon expansion 
around the upper critical dimension d c = 4. Our guess is that this new anomalous exponent 
is related to the renormalization of a composite operator not consider so far in the analysis 
of the Reggeon field theory, but this issue is beyond the scope of this paper, and will be 
studied elsewhere. 
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TABLES 



Dimension 


a 


«2 


Pw 


K 


z 


d= 1 


1.33(1) 


0.63(3) 


0.839(1) 


0.4336(4) 


1.58(1) 


d = 2 


0.97(1) 


0.385(5) 


0.550(5) 


0.33(1) 


1.765(10) 


d = 3 


0.51(1) 


0.09(2) 


0.27(1) 


0.22(1) 


1.90(5) 


d = 4 














2 



Summary of interface-growth critical exponents for the contact process obtained in our 
simulations (d = 1-3). Figures in parentheses denote uncertainties. 
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FIGURE CAPTIONS 

FIG. 1. Interface of the critical contact process in a system of 200 sites, shown at intervals 
of 5000 time units. 

FIG. 2. Scaled mean-square width versus reduced time for system sizes L = 500, 1000, 2000, 
and 5000. 

FIG. 3. Scaling plot of the height probability distribution (unnormalized) for L = 1000 at 
(from top to bottom on left-hand side) times 500, 1000, 2000, 5000, 10 4 , 2 x 10 4 , 5 x 10 4 , 
and 10 5 . 

FIG. 4. Scaling plot of the height probability distributions (unnormalized) for L = 1000 
(dashed lines) and L = 5000 (solid lines) at reduced times (from top to bottom on left-hand 
side) t = 0.142, 0.71, and 1.42. 

FIG. 5. Scaled height-height correlation function versus r/t 1 ^ for L = 1000. The topmost 
curve comprises collapsed data for times 500, 1000, 2000, 5000, and 10 4 ; below it lie results 
for t = 2 x 10 4 (solid curve), 5 x 10 4 (dotted curve), and 10 5 (dashed curve). 

FIG. 6. Growth of the mean-square gradient s(t) = (V/i) 2 in systems with L = 1000 (+), 
2000 (dashed line), and 5000 (solid line). 

FIG. 7. Scaled mean-square width versus reduced time in two dimensions; system sizes 
L = 32, 64, 128, and 256. 



FIG. 8. Growth of the mean-square gradient s(t) = (Vh) 2 in the two-dimensional system 
with L = 256. 

FIG. 9. Scaled height-height correlation function versus r/t 1 ^ for L = 256 in two dimensions. 
The curves (bottom to top) correspond to r = 2, 4, 8, and 16. 
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